# ==============================================================================
# DESCRIPTION ----
# ==============================================================================

# This script performs robustness checks for the event study analysis 

# ==============================================================================
# PREAMBLE, LIBRARIES, AND IMPORT ----
# ==============================================================================

#### Preamble ####
encoding = "UTF-8"
set.seed(999)

#### Libraries ####
library(tidyverse)
library(tidyquant)
library(eventstudies)

#### Import ####
# Import amakudari data
amakudari <- read_csv("data/amakudata.csv")

# ==============================================================================
# FILTER DATA TO VICE-MINISTERS ONLY ----
# ==============================================================================

# Filter data to vice-minister and assistant vice-minister retirements only
all <- amakudari %>% 
  filter(
    firm_type2_en == "Stock",
    position_level == "Vice-Minister" | 
      position_level == "Assistant Vice-Minister",
    #ministry != "Cabinet office",
    ministry != "Police",
    tse_code != -99,
    tse_code != 99) %>%
  select(date_ret, name, position_level, ministry_short, 
         firm_dest, firm_dest_en, post_dest, tse_code)

avm <- all %>% filter(position_level == "Assistant Vice-Minister")
vm <- all %>% filter(position_level == "Vice-Minister")

# ==============================================================================
# RETREIVE FINANCIAL DATA FROM YAHOO FINANCE ----
# ==============================================================================

# Create list of tickers
tickers <- unique(all$tse_code)
tickers <- paste(tickers, "T", sep = ".")

# Add some stock codes manually to increase sample of events
tickers <- append(tickers,'1812.T')
tickers <- append(tickers,'8411.T')

# Add Nikkei 225 Index and TOPIX for market model
tickers <- append(tickers, "^N225")

# Retrive financial data from Yahoo Finance using tidyquant
prices <- tq_get(tickers,
                 from = "2009-01-01",
                 to = "2019-12-31",
                 get = "stock.prices")

# Add percentage change column and transform long to wide
prices <- prices %>%
  group_by(symbol) %>% 
  mutate(pct_change = (adjusted/lag(adjusted) - 1) * 100) %>%
  ungroup(symbol) %>%
  filter(!is.na(pct_change)) %>%
  distinct(symbol, date, pct_change) %>%
  select(symbol, date, pct_change) %>%
  pivot_wider(names_from = symbol, values_from = pct_change)

# ==============================================================================
# ADD EVENT DATES FROM NEWS REPORTS TO ALL RETIREMENTS ----
# ==============================================================================

#### Import coded event dates ####
event_dates_avm <- read_csv("data/event_dates_avm.csv")
event_dates_vm <- read_csv("data/event_dates_vm.csv")
event_dates <- rbind(event_dates_avm, event_dates_vm) %>%
  rename(firm = firm_dest_english)

#### Create separate dataframes for each category #### 
# Create separate dataframe of outside directors
director <- event_dates %>% 
  filter((position == "internal director" | position == "outside director") & !is.na(when))

# Create separate dataframe of vice ministerial consultants
advisor_vm <- event_dates %>% 
  filter(position != "outside director" & !is.na(when)) %>%
  filter(position_level == "Vice-Minister" & !is.na(when))

# Create separate dataframe of all internal hires
internal <- event_dates %>% 
  filter((position == "advisor" | position == "executive" | position == "manager") & !is.na(when)) 

# Place in R date-time format: all
event_dates <- event_dates %>% 
  mutate(when = as.Date(when)) %>%
  filter(!is.na(when)) %>%
  select(ticker, when) %>%
  rename(name = ticker)

event_dates <- as.data.frame(event_dates)

# Place in R date-time format: outside directors
director_event_dates <- director %>% 
  mutate(when = as.Date(when)) %>%
  filter(!is.na(when)) %>%
  select(ticker, when) %>%
  rename(name = ticker)

director_event_dates <- as.data.frame(director_event_dates)

# Place in R date-time format: vice ministerial consultants
advisor_vm_event_dates <- advisor_vm %>% 
  mutate(when = as.Date(when)) %>%
  filter(!is.na(when)) %>%
  select(ticker, when) %>%
  rename(name = ticker)

advisor_vm_event_dates <- as.data.frame(advisor_vm_event_dates)

# Place in R date-time format: internal
internal_event_dates <- internal %>% 
  mutate(when = as.Date(when)) %>%
  filter(!is.na(when)) %>%
  select(ticker, when) %>%
  rename(name = ticker)

internal_event_dates <- as.data.frame(internal_event_dates)

# ==============================================================================
# PERFORM EVENT STUDIES USING CONSTANT MEAN RETURN MODEL ----
# ==============================================================================

# Convert stock data to zoo time series object
prices_ts <- read.zoo(prices)

# Define event window
event_window = 8

#### All appointments ####
cmr_all <- eventstudy(firm.returns = prices_ts,
                      event.list = event_dates,
                      event.window = event_window,
                      type = "None",
                      to.remap = TRUE,
                      remap = "cumsum",
                      inference = TRUE,
                      inference.strategy = "bootstrap")

plot(cmr_all)

#### Directors  ####
cmr_director <- eventstudy(firm.returns = prices_ts,
                          event.list = director_event_dates,
                          event.window = event_window,
                          type = "None",
                          to.remap = TRUE,
                          remap = "cumsum",
                          inference = TRUE,
                          inference.strategy = "bootstrap")

plot(cmr_director)

#### Internal hires ####
cmr_internal <- eventstudy(firm.returns = prices_ts,
                          event.list = internal_event_dates,
                          event.window = event_window,
                          type = "None",
                          to.remap = TRUE,
                          remap = "cumsum",
                          inference = TRUE,
                          inference.strategy = "bootstrap")

plot(cmr_internal)

#### Convert to tibble and clean data ####
mm_all <- as_tibble(cmr_all$result) %>%
  rename(mean = Mean, ci_low = `2.5%`, ci_high = `97.5%`) %>%
  mutate(event_day = row_number() - event_window) %>%
  mutate(id = "All")

mm_od <- as_tibble(cmr_director$result) %>%
  rename(mean = Mean, ci_low = `2.5%`, ci_high = `97.5%`) %>%
  mutate(event_day = row_number() - event_window) %>%
  mutate(id = "Directors")

mm_advisor <- as_tibble(cmr_internal$result) %>%
  rename(mean = Mean, ci_low = `2.5%`, ci_high = `97.5%`) %>%
  mutate(event_day = row_number() - event_window) %>%
  mutate(id = "Internal hires")

mm <- rbind(mm_all, mm_od, mm_advisor) %>%
  mutate(id = fct_relevel(id, "All", "Directors", "Internal hires"))

#### Plot cumulative abnormal returns and 95% confidence intervals ####
mm_plot = 
  ggplot(mm) +
  geom_line(aes(x = event_day, y = ci_low), linetype = "dashed") +
  geom_line(aes(x = event_day, y = ci_high), linetype = "dashed") +
  geom_ribbon(aes(x = event_day, ymin = ci_low, ymax = ci_high), 
              fill="steelblue2", alpha = 0.25) +
  geom_point(aes(x = event_day, y = mean)) +
  geom_line(aes(x = event_day, y = mean)) +
  geom_vline(xintercept = 0, linetype = "solid") +
  geom_hline(yintercept = 0, linetype = "solid") +
  xlab("Event day") + 
  ylab("Cumulative abnormal returns (%)") + 
  scale_x_continuous(limits = c(-8, 8), breaks=seq(-8, 8, 2)) +
  scale_y_continuous(limits = c(-7, 7), breaks=seq(-7, 7, 2)) +
  theme_classic() +
  theme(plot.title = element_text(hjust = 0.5)) +
  theme(axis.text=element_text(size = 10)) +
  theme(axis.text.x = element_text(size = 10)) +
  theme(legend.position = "none")

mm_plot = mm_plot + facet_wrap(vars(id), nrow = 3)

# FIGURE A30
ggsave("figures/a30_event_study_cmr.pdf", height = 6, width = 6)

# ==============================================================================
# CALCULATE STANDARD ERRORS USING T-TEST AND WILCOX RANK TEST ----
# ==============================================================================

#### All appointments: t-test ####
mm_all_t <- eventstudy(firm.returns = prices_ts,
                     event.list = event_dates,
                     event.window = event_window,
                     type = "marketModel",
                     to.remap = TRUE,
                     remap = "cumsum",
                     inference = TRUE,
                     inference.strategy = "classic",
                     model.args = list(market.returns = prices_ts$`^N225`))

#### Internal hires: t-test ####
mm_internal_t <- eventstudy(firm.returns = prices_ts,
                          event.list = internal_event_dates,
                          event.window = event_window,
                          type = "marketModel",
                          to.remap = TRUE,
                          remap = "cumsum",
                          inference = TRUE,
                          inference.strategy = "classic",
                          model.args = list(market.returns = prices_ts$`^N225`))

#### All appointments: wilcox ####
mm_all_wilcox <- eventstudy(firm.returns = prices_ts,
                     event.list = event_dates,
                     event.window = event_window,
                     type = "marketModel",
                     to.remap = TRUE,
                     remap = "cumsum",
                     inference = TRUE,
                     inference.strategy = "wilcox",
                     model.args = list(market.returns = prices_ts$`^N225`))

#### Internal hires: wilcox ####
mm_internal_wilcox <- eventstudy(firm.returns = prices_ts,
                          event.list = internal_event_dates,
                          event.window = event_window,
                          type = "marketModel",
                          to.remap = TRUE,
                          remap = "cumsum",
                          inference = TRUE,
                          inference.strategy = "wilcox",
                          model.args = list(market.returns = prices_ts$`^N225`))

#### Convert to tibble and clean data ####
mm_all_t <- as_tibble(mm_all_t$result) %>%
  rename(mean_median = Mean, ci_low = `2.5%`, ci_high = `97.5%`) %>%
  mutate(event_day = row_number() - event_window) %>%
  mutate(id = "All: T-test")

mm_internal_t <- as_tibble(mm_internal_t$result) %>%
  rename(mean_median = Mean, ci_low = `2.5%`, ci_high = `97.5%`) %>%
  mutate(event_day = row_number() - event_window) %>%
  mutate(id = "Internal hires: T-test")

mm_all_wilcox <- as_tibble(mm_all_wilcox$result) %>%
  rename(mean_median = Median, ci_low = `2.5%`, ci_high = `97.5%`) %>%
  mutate(event_day = row_number() - event_window) %>%
  mutate(id = "All: Wilcox")

mm_internal_wilcox <- as_tibble(mm_internal_wilcox$result) %>%
  rename(mean_median = Median, ci_low = `2.5%`, ci_high = `97.5%`) %>%
  mutate(event_day = row_number() - event_window) %>%
  mutate(id = "Internal hires: Wilcox")

mm <- rbind(mm_all_t, mm_internal_t, mm_all_wilcox, mm_internal_wilcox) %>%
  mutate(id = fct_relevel(id, "All: T-test", "Internal hires: T-test", 
                          "All: Wilcox", "Internal hires: Wilcox"))

#### Plot cumulate abnormal returns and 95% confidence intervals ####
mm_plot = 
  ggplot(mm) +
  geom_line(aes(x = event_day, y = ci_low), linetype = "dashed") +
  geom_line(aes(x = event_day, y = ci_high), linetype = "dashed") +
  geom_ribbon(aes(x = event_day, ymin = ci_low, ymax = ci_high), 
              fill="steelblue2", alpha = 0.25) +
  geom_point(aes(x = event_day, y = mean_median)) +
  geom_line(aes(x = event_day, y = mean_median)) +
  geom_vline(xintercept = 0, linetype = "solid") +
  geom_hline(yintercept = 0, linetype = "solid") +
  xlab("Event day") + 
  ylab("Cumulative abnormal returns (%)") + 
  scale_x_continuous(limits = c(-8, 8), breaks=seq(-8, 8, 2)) +
  scale_y_continuous(limits = c(-6, 6), breaks=seq(-6, 6, 2)) +
  theme_classic() +
  theme(plot.title = element_text(hjust = 0.5)) +
  theme(axis.text=element_text(size = 10)) +
  theme(axis.text.x = element_text(size = 10)) +
  theme(legend.position = "none")

mm_plot = mm_plot + facet_wrap(vars(id), nrow = 2)

#### Output ggplot ####
ggsave("figures/a31_event_study_classic_wilcox.pdf", height = 6, width = 6)

# ==============================================================================
# RECALCULATE RESULTS WITH OTHER EVENT WINDOWS ----
# ==============================================================================

#### Differing event windows ####
# Run event studies varying event windows
all_windows <- data.frame()
for (i in c(3, 4, 5, 6, 7, 9, 10, 15, 20)) {
  mm_advisor_windows <- eventstudy(firm.returns = prices_ts,
                                   event.list = internal_event_dates,
                                   event.window = i,
                                   type = "marketModel",
                                   to.remap = TRUE,
                                   remap = "cumsum",
                                   inference = TRUE,
                                   inference.strategy = "bootstrap",
                                   model.args = list(market.returns=prices_ts$`^N225`))
  
  mm_advisor_windows <- as_tibble(mm_advisor_windows$result) %>%
    rename(mean = Mean, ci_low = `2.5%`, ci_high = `97.5%`) %>%
    mutate(event_day = row_number() - i) %>%
    mutate(id = i)
  
  all_windows <- rbind(all_windows, mm_advisor_windows)
}

# Plot cumulate abnormal returns and 95% confidence intervals
mm_plot_windows = 
  ggplot(all_windows) +
  geom_line(aes(x = event_day, y = ci_low), linetype = "dashed") +
  geom_line(aes(x = event_day, y = ci_high), linetype = "dashed") +
  geom_ribbon(aes(x = event_day, ymin = ci_low, ymax = ci_high), 
              fill="steelblue2", alpha = 0.25) +
  geom_line(aes(x = event_day, y = mean)) +
  geom_vline(xintercept = 0, linetype = "solid") +
  geom_hline(yintercept = 0, linetype = "solid") +
  xlab("Event day") + 
  ylab("Cumulative abnormal returns (%)") + 
  scale_x_continuous(limits = c(-20, 20), breaks=seq(-20, 20, 5)) +
  scale_y_continuous(limits = c(-4, 8), breaks=seq(-4, 8, 2)) +
  theme_classic() +
  theme(plot.title = element_text(hjust = 0.5)) +
  theme(axis.text=element_text(size = 10)) +
  theme(axis.text.x = element_text(size = 10)) +
  theme(legend.position = "none")

mm_plot_windows = mm_plot_windows + facet_wrap(vars(id), nrow = 4)

#### Output ggplot ####
ggsave("figures/a32_event_study_event_windows.pdf", height = 6, width = 6)

# ==============================================================================
# TIME SHIFTED PLACEBO TEST ----
# ==============================================================================

#### Time-shifted placebo tests ####
all_placebos <- data.frame()
for (i in c(-200, -100, -50, -25, -10, -5, 0, 5, 10, 25, 50, 100, 200)) {
  
  placebo <- internal_event_dates %>% mutate(when = when + i)
  
  mm_advisor_placebos <- eventstudy(firm.returns = prices_ts,
                                    event.list = placebo,
                                    event.window = event_window,
                                    type = "marketModel",
                                    to.remap = TRUE,
                                    remap = "cumsum",
                                    inference = TRUE,
                                    inference.strategy = "bootstrap",
                                    model.args = list(market.returns=prices_ts$`^N225`))
  
  mm_advisor_placebos <- as_tibble(mm_advisor_placebos$result) %>%
    rename(mean = Mean, ci_low = `2.5%`, ci_high = `97.5%`) %>%
    mutate(event_day = row_number() - event_window) %>%
    mutate(id = i)
  
  all_placebos <- rbind(all_placebos, mm_advisor_placebos)
}

placebo_0 <- all_placebos %>% 
  filter(event_day == 0) %>%
  mutate(id = as.factor(id))

# Plot results of time-shifted placebo test
placebo_plot <- 
  ggplot(placebo_0, aes(x = id, y = mean)) +
  geom_hline(aes(yintercept = 0), linetype = 2, color = "grey") +
  geom_vline(aes(xintercept = "0"), linetype = 2, color = "grey") +
  geom_point(color = "steelblue2", size = 2) +
  geom_errorbar(aes(ymin = ci_low, ymax = ci_high), 
                color = "grey75", width = 0.1) +
  xlab("Number of days shifted from event date") +
  ylab("Mean event day abnormal return (%)") +
  scale_y_continuous(limits = c(-5, 5), 
                     breaks = scales::pretty_breaks(n = 10)) +
  theme_classic() +
  theme(plot.title = element_text(hjust = 0.5)) +
  theme(axis.text=element_text(size = 10)) +
  theme(axis.text.x = element_text(size = 10))

ggsave("figures/a29_event_study_placebo_test.pdf", height = 4, width = 6)

